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Abstract. Trapped degenerate dipolar Bose and Fermi gases of cylindrical symmetry 
r~~i with the polarization vector along the symmetry axis are only stable for the strength 

^ of dipolar interaction below a critical value. In the case of bosons, the stability of 

W) such a dipolar Bose-Einstcin condensate (BEC) is investigated for different strengths 

■y of contact and dipolar interactions using variational approximation and numerical 

^ solution of a mean-field model. In the disk shape, with the polarization vector 

^ perpendicular to the plane of the disk, the atoms experience an overall dipolar repulsion 

O^ and this fact should contribute to the stability. However, a complete numerical solution 

of the dynamics leads to the collapse of a strongly disk-shaped dipolar BEC due to 

the long-range anisotropic dipolar interaction. In the case of fermions, the stability 

of a trapped single-component degenerate dipolar Fermi gas is studied including the 

' ^ Hartree-Fock exchange and Brueckner-Goldstone correlation energies in the local- 

pi density approximation valid for a large number of atoms. Estimates for the maximum 

^ allowed number of polar Bose and Fermi molecules in BEC and degenerate Fermi gas 

are given. 
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1. Introduction 

After the experimental realization of Bose-Einstein condensate (BEC) of ^^Cr [U |2], 
^^^Dy ^, and ^^^Er [4j atoms and a single-component degenerate Fermi gas of ^^^Dy 
[3] atoms with large magnetic dipolar interaction, there has been new interest in the 
theoretical and experimental studies of degenerate gases. Polar Bose [6] and Fermi [7] 
molecules with much larger permanent electric dipole moment are also being considered 
for future experiments. The anisotropic long-range nonlocal dipolar interaction acting in 
all partial waves in these atoms and molecules is distinct from isotropic S'-wave contact 
interaction. Due to the anisotropic nonlocal nature of dipolar interaction, the stability 
of a dipolar BEC depends on the number of atoms, dipolar interaction and scattering 
length as well as, reasonably strongly and distinctly, on the trap geometry [21 El [ID]- 
For example, in the nondipolar case the spherically symmetric configuration is the most 
stable one |Tl] and accommodate the largest number of atoms, whereas in the dipolar 
case the disk-shaped configuration is the most stable one. 

One advantage in investigating the effect of dipolar interaction in a single- 
component degenerate dipolar Fermi gas over that in a dipolar BEC is the remarkable 
stability of the former. Due to three-body loss via molecule formation, a BEC can be 
easily destroyed. The three-body loss is highly suppressed in a degenerate Fermi gas due 
to Pauli repulsion among identical fermions and this system is unconditionally stable in 
the absence of dipolar interaction. Nevertheless, in the presence of dipolar interaction 
beyond a certain strength, both the trapped BEC [Bl El [ID] and the degenerate Fermi gas 
[121 [13] are unstable because of possible collapse. New mechanism and route to collapse 
open up in the presence of dipolar interaction because of atomic interaction in non-S* 
(angular momentum L ^ 0) waves. Dipolar BECs are immediately distinguishable 
from those with purely contact interaction by their strong and distinct shape and 
stability dependence on the trapping geometry. In the case of a trapped dipolar BEC 
of ^^Cr atoms, anisotropic D wave collapse has been studied both theoretically and 
experimentally [H]. Structure formation during the collapse has been studied in a 
single-component [lO] , in a binary dipolar BEC [I5] , and in a dipolar droplet bound in 
a nondipolar BEC [T^ . 

We undertake a systematic study of the stability of a trapped dipolar BEC and of a 
trapped degenerate dipolar Fermi gas, both with cylindrical symmetry. The polarization 
direction is always taken along the axial symmetry direction. Similar studies [TT] of 
nondipolar BECs illustrated the effect of the trapping symmetry and atomic contact 
attraction on stability. We present stability plots of a trapped dipolar BEC for varying 
number of atoms, scattering length, dipolar interaction, and trap symmetry. Such 
stability plots should aid in the future planning of experiments. The repulsive atomic 
contact interaction as well as strongly disk-shaped trap should facilitate stability and it 
is intuitively expected that these could lead to the absolute stability of the dipolar 
BEC and this possibility is supported by a Lagrangian variational analysis of the 
system. However, a complete dynamical treatment shows that for an increased dipolar 
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interaction, a cylindrically trapped dipolar BEC should collapse independent of the 
underlying trap symmetry and of a large atomic contact repulsion. In the case of polar 
molecules, under the usual condition of trapping, the number of molecules in a stable 
BEC could be very small, typically, less than a hundred. The number of molecules in 
a stable BEC can be increased by considering a significantly weaker trap with much 
reduced trapping frequency. 

We also study, using the real-time routine, the dynamical nature of local collapse in 
dipolar BECs induced by a change of dipolar interaction, as opposed to global collapse 
in nondipolar BEC induced by a change in contact interaction managed by the Feshbach 
resonance technique [T7] . The dipolar interaction can be controlled experimentally by a 
rotating orienting field [18] to initiate the collapse while maintaining a constant atomic 
interaction. The dipolar interaction is most prominent in disk and cigar shapes and we 
consider only these shapes. In both cases we evidence local collapse and not a global 
collapse to the trap center. In cigar shape, the dipolar BEC is first elongated and 
compressed on the axial symmetry direction and with further compression, centers of 
collapse develop along this axis and eventually the collapsed state with small pieces of 
BEC occupies a large extension along the symmetry axis. In disk shape, due to increased 
dipolar repulsion, the disk-shaped BEC first takes the shape of a donut and individual 
centers of collapse appear on this donut. Eventually, after the collapse and breaking up, 
small pieces of BEC occupy the whole volume of the entire initial disk. 

The stability of a trapped degenerate dipolar Fermi gas has been studied by several 
groups. Miyakawa et al. [ISj showed that the Fock exchange term causes a deformation 
of the Fermi surface in the presence of anisotropic dipolar interaction which plays a 
crucial role in the stability. They introduce a variational Wigner function to describe 
the deformation of the Fermi gas in phase space and use it to study the stability. Zhang 
and Yi [12] considered the deformed Fermi surface in the form of a ellipsoid and provided 
a variational and numerical solution of stability. Liu and Yin [19] further included a 
Brueckner-Goldstone (BG) correlation correction to the energy of a degenerate dipolar 
Fermi gas with deformed Fermi surface and studied its stability. However, we find 
contradiction among some of the studies and lack of physical plausibility of some of the 
results. In view of this, we find the necessity of a critical investigation of this problem. 



In Sec. 2J^ we present a brief description of the mean-field model of the dipolar BEC 
which we use in this work. We also present an analytical variational approximation of the 
same. The present theoretical formulation for the degenerate dipolar Fermi gas in local- 



density approximation (LDA) [20] is formulated in Sec. 2.2 An explicit contribution 
of Fock deformation to Fermi surface and of the BG correlation correction to energy 
and chemical potential is included in the LDA. We present a simple scaled form of the 



model which we use in numerical study. In Sec. 3.1, we present numerical results for 



the stability of a dipolar BEC and in Sec. 3.2 the same for the degenerate dipolar Fermi 



gas. Finally, in Sec. IV we present a summary of our study. 
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2. Analytical Consideration 

2.1. Mean-field model for dipolar BEC 

A dipolar BEC of N atoms, each of mass m satisfies the mean-field Gross-Pitaevskii 
(GP) equation [I] 



th^A 
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dt 
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with 0(r) the wave function, a the atomic scattering length, n(r) = 0^(r) the BEC 
density normalized as J n{v)dv = 1, A = oj^/cOp is the trap aspect ratio, to = {lo'^lOzY^^ 
is the average trap angular frequency, where Up and Uz are the angular frequencies of 
radial and axial traps. The dipolar interaction between two atoms at r and r' in ([I]) is 
taken as 



f/drifR) 



/io/irf (1 — 3cos^ 
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Kid(R), R = r-r' 
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where 6 is the angle between the vector R and the polarization direction z, /iq is the 
permeability of free space and fid is the magnetic dipole moment of each atom. In case 
of polar molecules, each of electric dipole moment d, the prefactor in ^ is (i^/(47reo), 
where eo is the permittivity of vacuum. To compare the strengths of atomic short-range 
and dipolar interactions, the dipolar interaction is often expressed in terms of the length 
scale add = ''^/Wo/irf/(127r/i^); for polar molecules add = md'^/{12Tieoh'^). Using this length 
scale, it is convenient to write the dipolar GP equation ([I]) in the following dimensionless 
form 
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In ([3]) energy, length, density r7,(r) and time t are expressed in units of oscillator energy 
hu, oscillator length Iq = Jh/mu, Iq^, and u'^, respectively. 

A partial understanding of the existence of a stable trapped BEC can be obtained 
from a variational approximation to (tsl) with the following ansatz 

1 
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where Wp and Wz are the variational widths along radial p and axial z directions. 
The effective energy per particle of a stationary state propagating as 0(r, t) ~ 
exp{—ifit/h)(f{r), where fi is the chemical potential, can be written as 



E = J dr[^\VMr)\' + ^(A-^V + X'/h')^'{r) 
+ 2naNip\r) + ^^ J drV'(r)v?^(r')Kid(r - r' 



(5) 
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With the density Q, the effective energy ([s]) becomes 

^^^Jf^ A^ N[a-a^^f{K)] 

2wl ^ Awl 2A2/3 ^ 4 ^ V^^wlw, ' ^ ' 

where k = Wp/wz and 

1 + 2^2 S/t^tanh-Vl - /^^ 
^('^^-y^^ (1 - «:2)3/2 • (7) 

The stable stationary states correspond to a global minimum of energy ^ for the 
ground state BEC. Maxima and saddle points have been identified to be excited 
stationary states [23]. The widths of the ground state BEC with Gaussian ansatz can 
be obtained by a minimization of energy (|6| by dE/dwz = dE/dwp = [21 1 



where 



Wp ^ 1 [g - add^(/t)/2] 2N 
A2/3 ^3+ ^2^ ^3^^' U 

^^^^,,,^1_ [a-a0^m 



wl v27r ^pW^z 
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The variational estimate for chemical potential yU is: 

^ 2m;2 ^ 4«;2 ^ 2A2/3 ^ 4 ^ ^2iwlwz ' ^ ^ 

2.2. Local density approximation for polarized dipolar fermions 

Based on a hydrodynamic description of the trapped degenerate dipolar gas of A^ single- 
component fermions each of mass m, we derived the following mean-field equation for 
this system 



^V^=[-'S^ + ^(^"V + AV3,.) 



im 2 



+ 



/iH + iV I f/dd(r-r')n(r')(ir' 



n(r), (13) 



/^H = ^[67r2iVn(r)]2/3, (14) 

and apphed it [26] to study statics and dynamics of a trapped degenerate ^^^Dy gas. 



In (13) /i is the chemical potential of the trapped gas, hh is the direct Hartree bulk 
chemical potential of a uniform gas with spherically symmetric Fermi surface, n{r) is 
the density and the dipolar interaction f/dd and the trapping parameters A and u are 
the same as in the case of bosons. 

The above consideration takes the Fermi surface to be spherical in the presence 
of the direct Hartree contribution to energy. However, in the presence of anisotropic 
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dipolar interaction the Fermi surface will be deformed p^ [T3l [T9l [27] and will take an 
ellipsoidal shape due to the exchange Fock interaction to energy. It has been shown 
that due to this deformation of the Fermi surface, the constant Fermi wave vector kp 
for the noninteracting Fermi gas gets deformed to [27] 

k{e) = kF + ^kU3cosH~l). (15) 

ovr 

This result is valid in the lowest order in the dipolar interaction strength add- To this 

order this deformation does not change the volume of the Fermi surface or the total 

energy or the bulk chemical potential. This leads to the following Hartree-Fock (HF) 

contribution to the bulk chemical potential in the second order of interaction strength 

[19l[27] 

f^HF = l^H- ^(67r2)i/3(3add)'{iVn(r)}^/l (16) 

In addition to the HF contribution to energy, there is another important 
contribution for the distorted Fermi surface. Kohn and Luttinger J2H] showed that, 
in the presence of an anisotropic Fermi surface, in addition to the HF contribution 
the bulk chemical potential another contribution from the Brueckner-Goldstone (BG) 
formula could be important. Lee and Yang [29] showed that the usual BG second-order 

(2) 

perturbative correction E^q to the ground-state energy of a dipolar Fermi gas has the 
following form: 

p(2) Y- (0|;7dd|m)(m|;7dd|0) .^„. 

^BG = 1^ ^ ^ > (17) 

where |0) and \m) are the ground and excited states of the noninteracting Fermi gas. 
This term has been evaluated by Monte Carlo integration and contributes the following 
to the bulk chemical potential [19] 

/igi = -y^(3add)^{iVr^(r)r/^ (18) 

so that the net bulk chemical potential including the HF (exchange) and BG 
(correlation) correction is given by 

fiHF+BG = f^HF - ^(3add)'{iV^(r)}^/l (19) 

9m 



The quantum pressure (gradient) term in (13) contributes much less than the 



dominant "Fermi energy" term Ti [67r^iVn(r)]^/^/(2?7i) and its neglect leads to the LDA 

[2ni[2ai2S] 

;. = ^ma;^(A-2/V + A^/V)+/io 

^3addiV^ r ^^ (20) 

m J 
where we have introduced a generic bulk chemical potential /iq for the uniform gas 
which can be hhi fJ^HF or fiHF+SG, respectively, for spherical Fermi surface, including 
HF exchange correction, or including HF exchange and BG correlation correction. 
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A convenient dimensionless form of the LDA (20) can be obtained in terms of the 



scahng length I = N^^^Iq where Iq = Jh/mco is the oscillator length. Then using 
the scaled quantities p = p/l,z = z/l,f = r/l,R = R/l,edd = '^N^^^io.dd/lo) = 



3N^/^add/l,fi = iJ,/{N^^^hu),n{r) = n{r)P, (20) can be written as 



;, = i(A-2/3p2 + x^'/'z') + ^[67r%(r)]2/3 



2' ^ '2' 

3 cos2 e 



+ e,,J j^, n{r')dr'-PelMr)f', (21) 

with Jdfn{r) = 1, where /3 = for Hartree, /3 = 28(67r2)i/3/i35 fo^ HF, and 



/3 = 28(67r2)V3/i35 + 14/9 for HF+BG. The scaled equation (21) is now independent 



of N and can be solved subject to the normalization condition to determine both the 



chemical potential and the density. This iV-independent universal scaling (21) of the 



hydrodynamical equation (13) has been possible only in the LDA (20) after neglecting 
the quantum pressure term. 

3. Numerical Result for Stability 

3.1. Dipolar BEC 

The dipolar BEC is known to be unstable to collapse [30] beyond a critical strength of 
dipolar interaction independent of trap symmetry even for repulsive contact interaction 
(positive values of scattering length a). In this study we consider repulsive contact 
interaction only. In variational equations (IS]) and ^ the contribution of the nonlocal 
dipolar interaction is reduced to a repulsive or attractive contribution contained in the 
terms —Odds' (k) ^t^^ —ciddh{K) similar to the contribution of the contact interaction 
through the scattering length a. In these equations the aspect ratio k =,< and > 1, 
corresponds to a spherical, cigar, and disk shape, respectively. Note that as n changes 
from (extreme cigar shape) to 00 (extreme disk shape), g{K,)/2 and h{K,) change 
from 1 to —2 while A changes from to 00. The variational approximation is an 
oversimplification of the actual state of affairs in case of dipolar interaction and fails to 
predict collapse in several cases. Improved variational ansatz [21] with time-dependent 
coupled Gaussian functions combined with spherical harmonics has been suggested and 
used in this context. In extreme cigar shape the variational contribution of the dipolar 
interaction is attractive and in extreme disk shape it is repulsive. In the disk shape one 
can argue erroneously that the system will always be stable for any number of atoms as 
parallel in-plane dipoles with dipole moments aligned perpendicular to the plane always 
repeal and should make the system more stable. The net effect of short-range atomic 
and long-range dipolar interactions contained in the square brackets of (tsl) and ^ has 
to be attractive (negative) beyond a critical value for the system to collapse. For a 
moderate to extreme disk-shaped trap k > 1 and the functions g^n) and h^n) in (Is]) 
and ([9]) are negative, and the contribution of add is positive (repulsive), thus prohibiting 
any collapse. For a moderate to extreme cigar-shaped trap k < 1 and the functions 
g^n) and /i(/t) are positive, and the net contribution of dipolar interaction in ([s]) and ^ 
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Figure 1. The critical value of the dipolar nonlinearity SNa^d versus the trap 
aspect ratio A = uJz/ujp for a stable dipolar BEC in the absence of atomic contact 
interaction (a = 0) for (a) positive and (b) negative values of the dipolar strength Odd 
from variational approximation (var) and numerical solution (num) of the mean-field 
equation. Length is expressed in units oi Iq = 1 /im. 



could be negative (attractive), thus leading to collapse beyond a critical value of dipolar 
interaction. In extreme cigar shape the functions g{K)/2 and /i(k) are positive but less 
than unity in magnitude. The quantities in square brackets in ([s]) and ^ are always 
positive (repulsive) for add < o,. Consequently, collapse in a cigar-shaped dipolar BEC 
in variational approximation is only possible for add > o, and in a strongly disk-shaped 
BEC it is never possible. Although, the dipolar BEC is unstable to collapse for large 
enough A^ for all values of parameters, the variational approximation can only predict 
collapse for add > cl- The chemical potential and sizes predicted by the variational 
approximation are in good agreement with the numerical results for sets of parameters 
where stable BECs exist. In addition, the variational approximation may predict stable 
states for parameters, for which the system actually collapses. 

We solve the GP equation ^ numerically by the spht-step Crank-Nicolson method 
[521 123] • The dipolar integral was evaluated in the Fourier momentum (k) space using 
convolution theorem as |33] 



dr'V,,ir - r')|0(r')P = ^-H-^M(k)^[|0p](k)}(r), (22) 

where J-'[] and J^^^{} are the Fourier transform (FT) and inverse FT, respectively. 
The FT of the dipole potential is analytically known [33]. The FT of density |0p is 
evaluated numerically by means of a standard fast FT (FFT) algorithm. The dipolar 
integral ^ is evaluated by the convolution rule (22). The inverse FT is taken by means 
of a standard FFT algorithm. The FFT is carried out in Cartesian coordinates and 
hence the GP equation is solved in three dimensions irrespective of the symmetry of 
the trapping potential. We use typically a space step of 0.1 and time step 0.001 and 
consider up to 512 points in space discretization in each direction. In our calculation 
we take the oscillator length /q = 1 A^m determined by the mass of an atom m and the 
angular frequency u. The stability of the BEC in the numerical solution of the GP 
equation is confirmed upon obtaining a converged result after time evolution for a long 
time: typically an interval of time At = 50 (units of u~^). If a BEC survives for such a 
long time, usually, it stays for ever and will be termed stable. 
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Figure 2. The critical number of atoms Nc versus tlie dipolar strength add for a stable 
dipolar BEC with atomic S'-wave scattering length (a) a = 150ao and (b) a — 50ao 
for trap aspect ratio A = 0.25 and 4 from variational approximation (v) and numerical 
solution (n) of the mean-field equation. Oscillator unit used in calculation is ^o = 1 
/im. 



First we consider the stability of the BEC after switching off the contact interaction 
by the Feshbach resonance technique (a = 0) [17j. This will show the effect of dipolar 
interaction alone. We consider both positive and negative values of the the dipolar 
strength add- The negative values of add can be achieved by a rotating orienting field 
[18] . When add is negative the dipolar interaction is attractive in disk shape and repulsive 
in cigar shape. So, for negative add, the BEC is more stable in cigar shape and less in 
disk shape. 

In figure [I] (a) and (b) we present the variational and numerical stability plots 
for a = and for positive and negative values of add, respectively, for different values 
of trap aspect ratio A. The numerical results show that for all values of A there can 
be collapse for both positive and negative add- The variational results fail to predict 
collapse for A > 5 in the case of positive add and for A < 0.15 in the case of negative 
add- In these cases, variational result predicts absolute stability of the dipolar BEC, 
whereas the numerical result predicts a critical value of number of atoms (or of the 
dipolar strength add) beyond which the dipolar BEC is unstable as shown in figures IT] 
(a) and (b). 

We also consider stability plots in the presence of a repulsive contact interaction 
(a > 0) for positive add- In the presence of one more variable in the system, the scattering 
length a, we present in figure |2] (a) the stability plot of the critical number of atoms 
Nc versus add for a = 150ao for A = 0.25 (cigar shape ) and A = 4 (disk shape), where 
ao is Bohr radius. The system collapses for A^ > Nc- The same plot for a = 50ao is 
shown in figure [2] (b). In both cases the variational result is also shown. The critical 
number of atoms Nc decreases as add increases for both disk and cigar shapes. In both 
cases variational result cannot predict collapse for A = 4 in the disk shape and predicts 
absolute stability for all N. For the cigar shape (A = 0.25) the variational result predicts 
collapse for add > ot and the variational stability line stays very close to the numerical 
stability line. 
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Figure 3. Dynamics to dipolar collapse: Isodensity contour of a disk-shaped dipolar 
BEC of 13000 ^^^Dy atoms, with scattering length a set to zero by the Feshbach 
resonance technique, at times (a) t = 0, (b) 2, (c) 3.5, and (d) 4. The trap parameters 
are w = 27r x 60 Hz, X — 8, l^ ~ 1 /^m, tg = 2.65 ms. The collapse is initiated with 
Odd = 16ao. At t = the dipolar strength a^d is jumped to add = 40ao by a rotating 
orienting field. The time is expressed in units of ip a-nd length in units of Iq. Density 
on surface of contour is 0.01 ij.m.~^. 



Wilson et al. [8] also considered stability plot of dipolar BEC using numerical and 
Gaussian variational methods. Their study is, at best, complimentary to the present 
one. The stability plot depends on three quantities for a fixed trap aspect ratio A: 
atomic scattering length a, dipolar strength a^d, and number of atoms A^. The study of 
Wilson et al, more appropriate for ^^Cr atoms, considers instability with a variation of 
scattering length a for fixed A^ and add- We consider instability against increase of A^ 
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for different add- So the present study should be useful to find the maximum number 
of atoms in a stable dipolar BEC for different dipolar atoms and should be useful in 
planning future experiments. We also show the domain of applicability of the Gaussian 
variational method: for moderate to extreme disk shape (A > 5) the Gaussian model 
always leads to stable dipolar BEC whereas the numerical analysis shows instability 
above a critical number of atoms. 




Figure 4. Dynamics to dipolar collapse: Isodensity contour of a cigar-shaped 
dipolar BEC of 1300 ^^''Dy atoms, with scattering length a set to zero by the Feshbach 
resonance technique, at times (a) t — 0, (b) 1, (c) 2, and (d) 4. The trap parameters 
are w = 27r x 60 Hz, A = 1/8, Iq = 1 /im, to = 2.65 ms. The collapse is initiated with 
Odd = 16ao. At t = the dipolar strength ajd is jumped to add = 48ao by a rotating 
orienting field. The time is expressed in units of to and length in units of Iq. Density 
on surface of contour is 0.005 fiin^^. 



Different experimental groups are trying to create a highly dipolar BEC of polar 
molecules [6]. The dipolar interaction strength of these molecules is typically Odd = 
2000ao [I] . The stability lines of figures [2] give an idea of the maximum number of 
these molecules in a trap with Iq = 1 fun. A preliminary calculation indicates these 
numbers to be only few tens. These numbers will increase if Iq is increased by reducing 
the trapping frequency u as the critical number of molecules for stability is inversely 
proportional to a/cJ. 

Now we study the nature and dynamics of collapse in case of a disk or cigar-shaped 
dipolar BEC with prominent dipolar interaction by real-time evolution. In the case 
of spherical symmetry, the net dipolar interaction tends to average out and in that 
case one can have a global collapse to the center of the dipolar BEC. To study a local 
collapse where a dipolar BEC collapses to different local centers we consider disk- and 
cigar-shaped dipolar BECs with a = and initiate the the collapse by jumping add to 
a larger value by a rotating orienting field |18j. In this effort, we consider ^^^Dy atoms 
with Odd = 130ao [3| and use the rotating orienting field to prepare an initial state 
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Figure 5. Chemical potential fi, energy E, and rms sizes (p), and (z), versus dipolar 
interaction strength add of a dipolar BEC of 10000 atoms from numerical solution 
and variational approximation of the mean-field GP equation for (a) scattering length 
a = 50ao and for trap aspect ratio A = 4 and for (b) a = 150ao and A = 0.25. The 
oscillator length is taken as /q = 1 /xm. The dark arrow on the right indicates the onset 

(c) 

of collapse in the numerical routine at Odd = cLm ■ 



with Odd = 16ao. The trap angular frequency u is taken as 27r x 60 Hz, so that the 
oscillator length Iq = 1 fim. The time scale to = 2.65 ms. We consider N = 13000 
atoms in a disk shaped trap with A = 8 and initiate the collapse by jumping add to 
40ao at time t = and study the subsequent real-time evolution of the dynamics. Four 
snapshots of the dipolar BEC at t/to = 0, 2, 3.5 and 4 are shown in figures [3] Due to 
in-plane dipolar repulsion, the atoms of the disk-shaped dipolar BEC of figure [3] (a) 
first moves away from center and takes a donut shape as in figure [3] (b) at t/to = 2. 
Then due to strong out-of-plane attraction the donut-shaped BEC starts the process of 
local collapse and breaks up into small pieces as shown in figure |3] (c) at t/to = 3.5. 
Finally, at t/to = 4 these small pieces undergo collapse and break into smaller pieces 
which eventually occupy the whole spatial extension of the disk as shown in figure [3] (d) 
at t/to = 4. 

Next we study the local collapse of a cigar-shaped dipolar BEC of A^ = 1300 ^^^Dy 
atoms in a trap with cj = 27r x 60 Hz, A = 1/8, to = 2.65 ms. The initial Odd is adjusted 
to 16ao by the rotating orienting field. A collapse in initiated by jumping add to 48ao at 
t = 0. The route to collapse is studied by real-time evolution of this dynamics as shown 
in figures |4j The dipolar BEC collapses on the symmetry axis, becomes elongated and 
eventually breaks up into small pieces at t/to = 4 along the symmetry axis. In both 
cases of disk and cigar shapes we have a local collapse instead of a global collapse to 
the global center of the dipolar BEC. 

The structures presented in figures |3] and |4] do not correspond to stationary states 
but to nonstationary nonequilibrium states. Although, the trapping potential has axial 
symmetry, we have not explicitly imposed this symmetry on the numerical solution of 
the GP equation. Consequently, the nonstationary nonequilibrium states of figures |3] 
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and |4] may exhibit a small violation of this symmetry. It was noted before by Wilson 
et al. [8j that the simple Gaussian wave function assumed in Q cannot approximate 
the local collapse shown, for example, in figures [3] and |4} Furthermore, Kreibich et al. 
[31] discussed this in detail and demonstrated that an improved ansatz for the stability 
fluctuations can overcome some of the restrictions of the Gaussian wave function. Also, 
the mean-field GP equation may not be fully adequate for the treatment of the multi- 
fragmented BEG, and a full many-body description of the dynamics might be necessary. 
We also study the chemical potential and root-mean-square (rms) sizes (p) and 
(z) of disk- and cigar-shaped dipolar BEGs for parameters corresponding to stability 
for < Odd < «dd where a^^ is the maximum value of a^^ beyond which the system 
collapses and no stationary state can be obtained. In the disk shape the variational 
approximation yields stable stationary state beyond this point, where the numerical 
routine predicts collapse. We take a = 50ao and lo = I /xm and show the present 
numerical and variational results for energy E, chemical potential fi, and rms sizes (z) 
and (/i) of 10000 atoms in figures [s] (a) and (b) for A = 4 (disk) and A = 0.25 (cigar), 
respectively. The agreement between numerical and variational results is satisfactory in 
all cases. The variational energies based on energy minimization are always larger than 
the numerical energies. 

3.2. Degenerate Fermi gas 



To find the stability plot of a single-component polarized Fermi gas we use (21) valid 



for a large number of fermions. We solve this equation numerically without further 



approximation. The dipolar nonlinearity in (20) could be very large for large A^, whereas 



the allowed values of the scaled dipolar interaction strength edd is a small number which 



makes the numerical treatment of the LDA (21) a relatively easy task. This inherent 
symmetry of the LDA of the dipolar Fermi gas can be turned to a good advantage in 
the calculation of the critical number of dipolar Fermi atoms in a trapped system. The 



solution of the scaled equation (21) with small nonlinearity (A^ = 1) can be related to 



the solution of (20) for large A^. We shall see that (21) permits solution for a maximum 
value ^crit of the dipolar strength ^dd- This can be used to find the actual critical number 
of atoms Nc in a trap for a specific dipolar interaction strength add using the scaling 
relation ^crit = 3A'c/^add/^o- 

In figure [6] (a) we plot the critical value of the scaled nonlinearity ^dd versus trap 
aspect ratio A for positive add including only the direct Hartree (H) contribution, direct 
Hartree plus exchange Fock (HF) contribution, and finally the HF plus correlation BG 
(HF-I-BG) contribution to the energy of fermions. The direct Hartree contribution 
to the energy is the largest, the inclusion of the exchange and correlation corrections 
reduce the energy. Hence, the Hartree contribution leads to the largest critical scaled 
nonlinearity, the HF and HF-I-BG contributions lead to smaller values of critical scaled 
nonlinearity as can be seen in figure M (a). In this figure we also plot the results of 
Miyakawa et al |T3] and of Zhang and Yi [12] . There is considerable discrepancy between 
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Figure 6. The critical value of the dipolar strength £crit = 3iV^/^add/^o versus the 
trap aspect ratio A for a normal dipolar Fermi gas for both (a) positive and (b) negative 

Add. 



these two previous calculations obtained using similar physical models, although the 
approximations and calculational procedure of these two studies are not identical. In 
view of this we consider the contribution of the different terms in the present model 
and compare with these previous calculations. The calculation of Zhang and Yi [12] 
including exchange is reasonably close to the present HF result up to A = 4 beyond 
which small discrepancy is noted. The calculation of Miyakawa et al [13] including the 
HF contribution to energy is quite distinct. The latter calculation including exchange 
agrees with the present direct Hartree result, and not with the present HF result, up 
to A = 2. The critical nonlinearity of Miyakawa et al including exchange for A > 2 is 
much larger than the present direct Hartree result, whereas after including exchange it 
should be smaller. The reason for this discrepancy is not fully clear, although it could 
be related to the approximate variational nature of this previous calculation [M]. For 
positive Odd the dipolar interaction is attractive in the cigar shape and repulsive in the 
disk shape and consequently, ^crit increases as A changes from 0.1 to 10. 

In figure [6] (b) we plot the same stabihty lines of figure [6] (a) for negative Odd- 
In this figure we also show the stability line of Zhang and Yi [12] including exchange. 
Considerable discrepancy is noted between the present HF model calculation and that 
of Ref. [12] for A < 2. For negative add the dipolar interaction is attractive in the 
disk shape and repulsive in the cigar shape and consequently, critical e decreases as A 
changes from 0.1 to 10. However, in the stability plot of Zhang and Yi the critical values 
of e for A = 0.1 and 10 are practically the same, which does not seem physically very 
plausible. In view of these crucial discrepancies, further studies are needed to unveil the 
truth. 

The stability lines of figure [6] (a) give an idea of the maximum number of polar Fermi 
molecules in a typical trap of oscillator length Iq = 1 fim.. Using the present HF+BG 
result, we can estimate the critical number of molecules for A = 0.1 and 10 using the 
critical values £crit = 1-1 and 2.5, respectively. Using a dipolar strength add = 2000ao 
\T], we obtain the critical number of molecules Nc ~ 1700 and 240,000, respectively, for 



Stability of trapped degenerate dipolar Bose and Fermi gases 15 

A = 0.1 and 10. As expected, these numbers are much larger than those in the case of 
bosonic polar molecules. The higher stability in case of the fermions is due to the Pauli 
exclusion principle among identical fermions. 

4. Summary and Conclusion 

Trapped nondipolar BECs with repulsive atomic interaction as well as trapped 
degenerate single-component Fermi gas are unconditionally stable. The BECs are 
subject to collapse instability for atomic attraction above a critical value [H]. Trapped 
dipolar BECs with a sufficiently strong dipolar interaction are more fragile and could be 
unstable to collapse independent of the underlying trap symmetry even for repulsive 
contact interaction [HI [9]. The condition of stability of the dipolar BECs is fully 
distinct. We studied the conditions of stability of dipolar BECs in axially-symmetric 
trap for varying trap symmetry and varying strengths of contact and dipolar atomic 
interactions. In all cases collapse is possible for a sufficiently large dipolar interaction. 
We demonstrate the local nature of dipolar collapse by real-time simulation of the 
collapse dynamics, where centers of local collapse develop all over the condensate and 
small fragments of the condensate collapse to these local centers, as opposed to a global 
center. We estimate the maximum number of bosonic polar molecules in a stable 
condensate. For a dipolar interaction strength a^d = 2000ao, the maximum number 
of molecules in a trap of oscillator length Iq = 1 fun is typically few tens and this 
number can be increased with a weaker trap. 

We also studied the stability of a degenerate Fermi gas in stability plots. 
Independent of the underlying trap symmetry, this system will collapse for a sufficiently 
strong dipolar interaction [13j. The maximum number of Fermi atoms allowed in a 
stable system is usually much larger than that permitted in a dipolar BEC under similar 
conditions. We also estimate the maximum number of fermionic polar molecules in a 
stable system. For a typical dipolar interaction strength add = 2000ao, the maximum 
number of fermionic molecules in a trap of oscillator length /q = 1 /^ni is larger than the 
maximum number of bosonic molecules by several orders of magnitude. 
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